{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "LOWANK49Pi27"
   },
   "source": [
    "# 第十七章 潜在语义分析"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "1.单词向量空间模型通过单词的向量表示文本的语义内容。以单词-文本矩阵$X$为输入，其中每一行对应一个单词，每一列对应一个文本，每一个元素表示单词在文本中的频数或权值（如TF-IDF）\n",
    "$$X = \\left[ \\begin{array} { c c c c } { x _ { 11 } } & { x _ { 12 } } & { \\cdots } & { x _ { 1 n } } \\\\ { x _ { 21 } } & { x _ { 22 } } & { \\cdots } & { x _ { 2 n } } \\\\ { \\vdots } & { \\vdots } & { } & { \\vdots } \\\\ { x _ { m 1 } } & { x _ { m 2 } } & { \\cdots } & { x _ { m n } } \\end{array} \\right]$$\n",
    "单词向量空间模型认为，这个矩阵的每一列向量是单词向量，表示一个文本，两个单词向量的内积或标准化内积表示文本之间的语义相似度。\n",
    "\n",
    "2.话题向量空间模型通过话题的向量表示文本的语义内容。假设有话题文本矩阵$$Y = \\left[ \\begin{array} { c c c c } { y _ { 11 } } & { y _ { 12 } } & { \\cdots } & { y _ { 1 n } } \\\\ { y _ { 21 } } & { y _ { 22 } } & { \\cdots } & { y _ { 2 n } } \\\\ { \\vdots } & { \\vdots } & { } & { \\vdots } \\\\ { y _ { k 1 } } & { y _ { k 2 } } & { \\cdots } & { y _ { k n } } \\end{array} \\right]$$\n",
    "其中每一行对应一个话题，每一列对应一个文本每一个元素表示话题在文本中的权值。话题向量空间模型认为，这个矩阵的每一列向量是话题向量，表示一个文本，两个话题向量的内积或标准化内积表示文本之间的语义相似度。假设有单词话题矩阵$T$\n",
    "$$T = \\left[ \\begin{array} { c c c c } { t _ { 11 } } & { t _ { 12 } } & { \\cdots } & { t _ { 1 k } } \\\\ { t _ { 21 } } & { t _ { 22 } } & { \\cdots } & { t _ { 2 k } } \\\\ { \\vdots } & { \\vdots } & { } & { \\vdots } \\\\ { t _ { m 1 } } & { t _ { m 2 } } & { \\cdots } & { t _ { m k } } \\end{array} \\right]$$ \n",
    "其中每一行对应一个单词，每一列对应一个话题，每一个元素表示单词在话题中的权值。\n",
    "\n",
    "给定一个单词文本矩阵$X$\n",
    "$$X = \\left[ \\begin{array} { c c c c } { x _ { 11 } } & { x _ { 12 } } & { \\cdots } & { x _ { 1 n } } \\\\ { x _ { 21 } } & { x _ { 22 } } & { \\cdots } & { x _ { 2 n } } \\\\ { \\vdots } & { \\vdots } & { } & { \\vdots } \\\\ { x _ { m 1 } } & { x _ { m 2 } } & { \\cdots } & { x _ { m n } } \\end{array} \\right]$$\n",
    "\n",
    "潜在语义分析的目标是，找到合适的单词-话题矩阵$T$与话题文本矩阵$Y$，将单词文本矩阵$X$近似的表示为$T$与$Y$的乘积形式。\n",
    "$$X \\approx T Y$$\n",
    "\n",
    "等价地，潜在语义分析将文本在单词向量空间的表示X通过线性变换$T$转换为话题向量空间中的表示$Y$。\n",
    "\n",
    "潜在语义分析的关键是对单词-文本矩阵进行以上的矩阵因子分解（话题分析）\n",
    "\n",
    "3.潜在语义分析的算法是奇异值分解。通过对单词文本矩阵进行截断奇异值分解，得到\n",
    "$$X \\approx U _ { k } \\Sigma _ { k } V _ { k } ^ { T } = U _ { k } ( \\Sigma _ { k } V _ { k } ^ { T } )$$\n",
    "\n",
    "矩阵$U_k$表示话题空间，矩阵$( \\Sigma _ { k } V _ { k } ^ { T } )$是文本在话题空间的表示。\n",
    "\n",
    "4.非负矩阵分解也可以用于话题分析。非负矩阵分解将非负的单词文本矩阵近似分解成两个非负矩阵$W$和$H$的乘积，得到\n",
    "$$X \\approx W H$$\n",
    "\n",
    "矩阵$W$表示话题空间，矩阵$H$是文本在话题空间的表示。\n",
    "\n",
    "非负矩阵分解可以表为以下的最优化问题：\n",
    "$$\\left. \\begin{array} { l } { \\operatorname { min } _ { W , H } \\| X - W H \\| ^ { 2 } } \\\\ { \\text { s.t. } W , H \\geq 0 } \\end{array} \\right.$$\n",
    "非负矩阵分解的算法是迭代算法。乘法更新规则的迭代算法，交替地对$W$和$H$进行更新。本质是梯度下降法，通过定义特殊的步长和非负的初始值，保证迭代过程及结果的矩阵$W$和$H$均为非负。\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "0WD7XVWRPkX1"
   },
   "source": [
    "-----\n",
    "**LSA** 是一种无监督学习方法，主要用于文本的话题分析，其特点是通过矩阵分解发现文本与单词之间的基于话题的语义关系。也称为潜在语义索引（Latent semantic indexing, LSI）。\n",
    "\n",
    "LSA 使用的是非概率的话题分析模型。将文本集合表示为**单词-文本矩阵**，对单词-文本矩阵进行**奇异值分解**，从而得到话题向量空间，以及文本在话题向量空间的表示。\n",
    "\n",
    "**非负矩阵分解**（non-negative matrix factorization, NMF）是另一种矩阵的因子分解方法，其特点是分解的矩阵非负。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "P1sWKgTGQ7r-"
   },
   "source": [
    "## 单词向量空间  \n",
    "word vector space model"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "CqXj1777RM8y"
   },
   "source": [
    "给定一个文本，用一个向量表示该文本的”语义“， 向量的**每一维对应一个单词**，其数值为该单词在该文本中出现的频数或权值；基本假设是文本中所有单词的出现情况表示了文本的语义内容，文本集合中的每个文本都表示为一个向量，存在于一个向量空间；向量空间的度量，如内积或标准化**内积**表示文本之间的**相似度**。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "3HVCXf6CSmTT"
   },
   "source": [
    "给定一个含有$n$个文本的集合$D=({d_{1}, d_{2},...,d_{n}})$，以及在所有文本中出现的$m$个单词的集合$W=({w_{1},w_{2},...,w_{m}})$. 将单词在文本的出现的数据用一个单词-文本矩阵（word-document matrix）表示，记作$X$:\n",
    "\n",
    "$\n",
    "X = \\begin{bmatrix}\n",
    "x_{11} &  x_{12}&  x_{1n}& \\\\ \n",
    "x_{21}&  x_{22}&  x_{2n}& \\\\ \n",
    "\\vdots &  \\vdots &  \\vdots & \\\\ \n",
    "x_{m1}&  x_{m2}&  x_{mn}& \n",
    "\\end{bmatrix}\n",
    "$\n",
    "\n",
    "这是一个$m*n$矩阵，元素$x_{ij}$表示单词$w_{i}$在文本$d_{j}$中出现的频数或权值。由于单词的种类很多，而每个文本中出现单词的种类通常较少，所有单词-文本矩阵是一个稀疏矩阵。\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "K2ncB3cde1Ab"
   },
   "source": [
    "权值通常用单词**频率-逆文本率**（term frequency-inverse document frequency, TF-IDF）表示：\n",
    "\n",
    "$TF-IDF(t, d ) = TF(t, d) * IDF(t)$,  \n",
    "\n",
    "其中，$TF(t,d)$为单词$t$在文本$d$中出现的概率，$IDF(t)$是逆文本率，用来衡量单词$t$对表示语义所起的重要性， \n",
    "\n",
    "$IDF(t) = log(\\frac{len(D)}{len(t \\in D) + 1})$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "bpu7MycIgu65"
   },
   "source": [
    "单词向量空间模型的优点是**模型简单，计算效率高**。因为单词向量通常是稀疏的，单词向量空间模型也有一定的局限性，体现在内积相似度未必能够准确表达两个文本的语义相似度上。因为自然语言的单词具有一词多义性（polysemy）及多词一义性（synonymy）。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "ns5wncZohn-z"
   },
   "source": [
    "## 话题向量空间"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "mmZpPHIdhrAy"
   },
   "source": [
    "**1. 话题向量空间**：\n",
    "\n",
    "给定一个含有$n$个文本的集合$D=({d_{1}, d_{2},...,d_{n}})$，以及在所有文本中出现的$m$个单词的集合$W=({w_{1},w_{2},...,w_{m}})$. 可以获得其单词-文本矩阵$X$:  \n",
    "\n",
    "$\n",
    "X = \\begin{bmatrix}\n",
    "x_{11} &  x_{12}&  x_{1n}& \\\\ \n",
    "x_{21}&  x_{22}&  x_{2n}& \\\\ \n",
    "\\vdots &  \\vdots &  \\vdots & \\\\ \n",
    "x_{m1}&  x_{m2}&  x_{mn}& \n",
    "\\end{bmatrix}\n",
    "$\n",
    "\n",
    "\n",
    "假设所有文本共含有$k$个话题。假设每个话题由一个定义在单词集合$W$上的$m$维向量表示，称为话题向量，即：  \n",
    "$t_{l} = \\begin{bmatrix}\n",
    "t_{1l}\\\\ \n",
    "t_{2l}\\\\ \n",
    "\\vdots \\\\ \n",
    "t_{ml}\\end{bmatrix},  l=1,2,...,k$\n",
    "\n",
    "其中$t_{il}$单词$w_{i}$在话题$t_{l}$的权值，$i=1,2,...,m$, 权值越大，该单词在该话题中的重要程度就越高。这$k$个话题向量 $t_{1},t_{2},...,t_{k}$张成一个话题向量空间（topic vector space）， 维数为$k$.**话题向量空间是单词向量空间的一个子空间**。\n",
    "\n",
    "话题向量空间$T$：  \n",
    "\n",
    "\n",
    "$\n",
    "T = \\begin{bmatrix}\n",
    "t_{11} &  t_{12}&  t_{1k}& \\\\ \n",
    "t_{21}&  t_{22}&  t_{2k}& \\\\ \n",
    "\\vdots &  \\vdots &  \\vdots & \\\\ \n",
    "t_{m1}&  t_{m2}&  t_{mk}& \n",
    "\\end{bmatrix}\n",
    "$  \n",
    "\n",
    "矩阵$T$,称为**单词-话题矩阵**。 $T = [t_{1},  t_{2}, ..., t_{k}]$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "Oc1c3JcKlTBD"
   },
   "source": [
    "**2. 文本在话题向量空间中的表示**  ：\n",
    "\n",
    "考虑文本集合$D$的文本$d_{j}$, 在单词向量空间中由一个向量$x_{j}$表示，将$x_{j}$投影到话题向量空间$T$中，得到话题向量空间的一个向量$y_{j}$, $y_{j}$是一个$k$维向量：  \n",
    "\n",
    "$y_{j} = \\begin{bmatrix}\n",
    "y_{1j}\\\\ \n",
    "y_{2j}\\\\ \n",
    "\\vdots \\\\ \n",
    "y_{kj}\\end{bmatrix},  j=1,2,...,n$ \n",
    "\n",
    "其中，$y_{lj}$是文本$d_{j}$在话题$t_{l}$中的权值， $l = 1,2,..., k$, 权值越大，该话题在该文本中的重要程度就越高。  \n",
    "\n",
    "矩阵$Y$ 表示话题在文本中出现的情况，称为话题-文本矩阵（topic-document matrix）,记作：  \n",
    "\n",
    "$\n",
    "Y = \\begin{bmatrix}\n",
    "y_{11} &  y_{12}&  y_{1n}& \\\\ \n",
    "y_{21}&  y_{22}&  y_{2n}& \\\\ \n",
    "\\vdots &  \\vdots &  \\vdots & \\\\ \n",
    "y_{k1}&  y_{k2}&  y_{kn}& \n",
    "\\end{bmatrix}\n",
    "$  \n",
    "\n",
    "也可写成： $Y = [y_{1}, y_{2} ..., y_{n}]$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "YcU3xwYindTo"
   },
   "source": [
    "**3. 从单词向量空间到话题向量空间的线性变换**：  \n",
    "\n",
    "如此，单词向量空间的文本向量$x_{j}$可以通过他在话题空间中的向量$y_{j}$近似表示，具体地由$k$个话题向量以$y_{j}$为系数的线性组合近似表示：  \n",
    "\n",
    "$x_{j} = y_{1j}t_{1} + y_{2j}t_{2} + ... + y_{yj}t_{k}, j = 1,2,..., n$  \n",
    "\n",
    "所以，单词-文本矩阵$X$可以近似的表示为单词-话题矩阵$T$与话题-文本矩阵$Y$的乘积形式。\n",
    "\n",
    "$X \\approx TY$  \n",
    "\n",
    "直观上，潜在语义分析是将单词向量空间的表示通过线性变换转换为在话题向量空间中的表示。这个线性变换由矩阵因子分解式的形式体现。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "Cu4JekfXFMqs"
   },
   "source": [
    "### 潜在语义分析算法  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "0awNXCy1Gw0K"
   },
   "source": [
    "潜在语义分析利用矩阵奇异值分解，具体地，对单词-文本矩阵进行奇异值分解，将其左矩阵作为话题向量空间，将其对角矩阵与右矩阵的乘积作为文本在话题向量空间的表示。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "otq3HMu5HVoK"
   },
   "source": [
    "给定一个含有$n$个文本的集合$D=({d_{1}, d_{2},...,d_{n}})$，以及在所有文本中出现的$m$个单词的集合$W=({w_{1},w_{2},...,w_{m}})$. 可以获得其单词-文本矩阵$X$:  \n",
    "$\n",
    "X = \\begin{bmatrix}\n",
    "x_{11} &  x_{12}&  x_{1n}& \\\\ \n",
    "x_{21}&  x_{22}&  x_{2n}& \\\\ \n",
    "\\vdots &  \\vdots &  \\vdots & \\\\ \n",
    "x_{m1}&  x_{m2}&  x_{mn}& \n",
    "\\end{bmatrix}\n",
    "$\n",
    "\n",
    "\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "mwNGRDgrHmmV"
   },
   "source": [
    "**截断奇异值分解**：\n",
    "\n",
    "潜在语义分析根据确定的话题数$k$对单词-文本矩阵$X$进行截断奇异值分解：  \n",
    "\n",
    "$\n",
    "X \\approx U_{k}\\Sigma _{k}V_{k}^{T} = \\begin{bmatrix}\n",
    "\\mu _{1} &  \\mu _{2}&  \\cdots & \\mu _{k}\n",
    "\\end{bmatrix}\\begin{bmatrix}\n",
    "\\sigma_{1} &  0&  0& 0\\\\ \n",
    " 0&  \\sigma_{2}&  0& 0\\\\ \n",
    " 0&  0&  \\ddots & 0\\\\ \n",
    " 0&  0&  0& \\sigma_{k}\n",
    "\\end{bmatrix}\\begin{bmatrix}\n",
    "v_{1}^{T}\\\\ \n",
    "v_{2}^{T}\\\\ \n",
    "\\vdots \\\\ \n",
    "v_{k}^{T}\\end{bmatrix}\n",
    "$\n",
    "\n",
    "矩阵$U_{k}$的每一个列向量 $u_{1}, u_{2},..., u_{k}$ 表示一个话题，称为**话题向量**。由这 $k$ 个话题向量张成一个子空间：  \n",
    "\n",
    "$\n",
    "U_{k} = \\begin{bmatrix}\n",
    "u_{1} &  u_{2}&  \\cdots & u_{k}\n",
    "\\end{bmatrix}\n",
    "$\n",
    "\n",
    "称为**话题向量空间**。  \n",
    "\n",
    "综上， 可以通过对单词-文本矩阵的奇异值分解进行潜在语义分析：  \n",
    "\n",
    "$ X \\approx U_{k} \\Sigma_{k} V_{k}^{T} = U_{k}(\\Sigma_{k}V_{k}^{T})$  \n",
    "\n",
    "得到话题空间 $U_{k}$ ， 以及文本在话题空间的表示($\\Sigma_{k}V_{k}^{T}$). "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "UTNHyq8mK8l5"
   },
   "source": [
    "### 非负矩阵分解算法"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "RQvqaMDYK_jf"
   },
   "source": [
    "非负矩阵分解也可以用于话题分析。对单词-文本矩阵进行非负矩阵分解，将**其左矩阵作为话题向量空间**，将其**右矩阵作为文本在话题向量空间的表示**。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "ApM8tE3MLqpP"
   },
   "source": [
    "#### 非负矩阵分解  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "glMwmkiwLyIn"
   },
   "source": [
    "若一个矩阵的索引元素非负，则该矩阵为非负矩阵。若$X$是非负矩阵，则： $X >= 0$.  \n",
    "\n",
    "给定一个非负矩阵$X$, 找到两个非负矩阵$W >= 0$ 和 $H>= 0$, 使得：  \n",
    "\n",
    "$ X \\approx WH$\n",
    "\n",
    "即非负矩阵$X$分解为两个非负矩阵$W$和$H$的乘积形式，成为非负矩阵分解。因为$WH$与$X$完全相等很难实现，所以只要求近似相等。  \n",
    "\n",
    "假设非负矩阵$X$是$m\\times n$矩阵，非负矩阵$W$和$H$分别为 $m\\times k$ 矩阵和 $k\\times n$ 矩阵。假设 $k < min(m, n)$ 即$W$ 和 $H$ 小于原矩阵 $X$, 所以非负矩阵分解是对原数据的压缩。\n",
    "\n",
    "称 $W$ 为基矩阵， $H$ 为系数矩阵。非负矩阵分解旨在用较少的基向量，系数向量来表示为较大的数据矩阵。\n",
    "\n",
    "令 $W = \\begin{bmatrix}\n",
    "w_{1} &  w_{2}&  \\cdots& w_{k} \n",
    "\\end{bmatrix}$\n",
    "为话题向量空间， $w_{1}, w_{2}, ..., w_{k}$ 表示文本集合的 $k$ 个话题， 令 $H = \\begin{bmatrix}\n",
    "h_{1} &  h_{2}&  \\cdots& h_{n} \n",
    "\\end{bmatrix}$\n",
    "为文本在话题向量空间的表示， $h_{1}, h_{2},..., h_{n}$ 表示文本集合的 $n$ 个文本。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "1DcvVSR0N_CF"
   },
   "source": [
    "##### 算法"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "hvZyHT85O5qt"
   },
   "source": [
    "非负矩阵分解可以形式化为最优化问题求解。可以利用平方损失或散度来作为损失函数。\n",
    "\n",
    "目标函数 $|| X - WH ||^{2}$ 关于 $W$ 和 $H$ 的最小化，满足约束条件 $W, H >= 0$, 即：  \n",
    "\n",
    "$\\underset{W,H}{min} || X - WH ||^{2}$  \n",
    "\n",
    "\n",
    "$s.t.  W, H >= 0$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "-zIGS1AEQdWp"
   },
   "source": [
    "乘法更新规则： \n",
    "\n",
    "\n",
    "$W_{il} \\leftarrow W_{il}\\frac{(XH^{T})_{il}}{(WHH^{T})_{il}}$  （17.33）\n",
    "\n",
    "\n",
    "$H_{lj} \\leftarrow H_{lj}\\frac{(W^{T}X)_{lj}}{(W^{T}WH)_{lj}}$  （17.34）\n",
    "\n",
    "\n",
    "选择初始矩阵 $W$ 和 $H$ 为非负矩阵，可以保证迭代过程及结果的矩阵 $W$ 和 $H$ 非负。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "MeiA0REkRpRi"
   },
   "source": [
    "**算法 17.1 （非负矩阵分解的迭代算法）**\n",
    "\n",
    "输入： 单词-文本矩阵 $X >= 0$, 文本集合的话题个数 $k$, 最大迭代次数 $t$；  \n",
    "输出： 话题矩阵 $W$, 文本表示矩阵 $H$。  \n",
    "\n",
    "**1)**. 初始化\n",
    "\n",
    "$W>=0$, 并对 $W$ 的每一列数据归一化；  \n",
    "$H>=0$；\n",
    "\n",
    "**2)**. 迭代  \n",
    "\n",
    "对迭代次数由1到$t$执行下列步骤：  \n",
    "a. 更新$W$的元素，对 $l$ 从1到 $k,i$从1到$m$按(17.33)更新 $W_{il}$;   \n",
    "a. 更新$H$的元素，对 $l$ 从1到 $k,j$从1到$m$按(17.34)更新 $H_{lj}$; "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "rIw6a0HITg08"
   },
   "source": [
    "### 图例 17.1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "colab": {},
    "colab_type": "code",
    "id": "0hPH9VEMPVGu"
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "from sklearn.decomposition import TruncatedSVD"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 125
    },
    "colab_type": "code",
    "id": "kjHirYzQWItl",
    "outputId": "0e6f2615-6a0b-4c4f-e74c-559727519eab"
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[2, 0, 0, 0],\n",
       "       [0, 2, 0, 0],\n",
       "       [0, 0, 1, 0],\n",
       "       [0, 0, 2, 3],\n",
       "       [0, 0, 0, 1],\n",
       "       [1, 2, 2, 1]])"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "X = [[2, 0, 0, 0], [0, 2, 0, 0], [0, 0, 1, 0], [0, 0, 2, 3], [0, 0, 0, 1], [1, 2, 2, 1]]\n",
    "X = np.asarray(X);X"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "colab": {},
    "colab_type": "code",
    "id": "I2yFnNJKWcPP"
   },
   "outputs": [],
   "source": [
    "# 奇异值分解\n",
    "U,sigma,VT=np.linalg.svd(X)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 233
    },
    "colab_type": "code",
    "id": "ollDH_QNXAdY",
    "outputId": "8d599d82-4bae-4047-941c-8ca524142349"
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[-7.84368672e-02, -2.84423033e-01,  8.94427191e-01,\n",
       "        -2.15138396e-01, -6.45002451e-02, -2.50012770e-01],\n",
       "       [-1.56873734e-01, -5.68846066e-01, -4.47213595e-01,\n",
       "        -4.30276793e-01, -1.29000490e-01, -5.00025540e-01],\n",
       "       [-1.42622354e-01,  1.37930417e-02, -1.25029761e-16,\n",
       "         6.53519444e-01,  3.88575115e-01, -6.33553733e-01],\n",
       "       [-7.28804669e-01,  5.53499910e-01, -2.24565656e-16,\n",
       "        -1.56161345e-01, -3.23288048e-01, -1.83248673e-01],\n",
       "       [-1.47853320e-01,  1.75304609e-01,  8.49795536e-18,\n",
       "        -4.87733411e-01,  8.40863653e-01,  4.97204799e-02],\n",
       "       [-6.29190197e-01, -5.08166890e-01, -1.60733896e-16,\n",
       "         2.81459486e-01,  1.29000490e-01,  5.00025540e-01]])"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "U"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 35
    },
    "colab_type": "code",
    "id": "lmxB_JViXFAF",
    "outputId": "9cdcba5c-74f8-42e4-9a4b-fad1c06b83cd"
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([4.47696617, 2.7519661 , 2.        , 1.17620428])"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sigma"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 161
    },
    "colab_type": "code",
    "id": "AiXURUScXMsj",
    "outputId": "d6fc0576-8bf1-491e-caed-02035079f0d3"
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[-1.75579600e-01, -3.51159201e-01, -6.38515454e-01,\n",
       "        -6.61934313e-01],\n",
       "       [-3.91361272e-01, -7.82722545e-01,  3.79579831e-02,\n",
       "         4.82432341e-01],\n",
       "       [ 8.94427191e-01, -4.47213595e-01, -2.23998094e-16,\n",
       "         5.45344065e-17],\n",
       "       [-1.26523351e-01, -2.53046702e-01,  7.68672366e-01,\n",
       "        -5.73674125e-01]])"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "VT"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 53
    },
    "colab_type": "code",
    "id": "DKOxld5lXRCK",
    "outputId": "0832796a-9952-4f39-e1ef-79347c495d16"
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "TruncatedSVD(algorithm='randomized', n_components=3, n_iter=7, random_state=42,\n",
       "             tol=0.0)"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# 截断奇异值分解\n",
    "\n",
    "svd = TruncatedSVD(n_components=3, n_iter=7, random_state=42)\n",
    "svd.fit(X)  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 35
    },
    "colab_type": "code",
    "id": "btnGrF0LXzZI",
    "outputId": "ba85127a-4fa6-4092-828c-9036c47f82f6"
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[0.39945801 0.34585056 0.18861789]\n"
     ]
    }
   ],
   "source": [
    "print(svd.explained_variance_ratio_)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 35
    },
    "colab_type": "code",
    "id": "F1hSe5NxX1zw",
    "outputId": "b0d0b87d-195b-4653-a857-48ff1eca887a"
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.933926460028446\n"
     ]
    }
   ],
   "source": [
    "print(svd.explained_variance_ratio_.sum())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 35
    },
    "colab_type": "code",
    "id": "cV4L2i9WX30R",
    "outputId": "6c313215-d095-41b2-a384-3dc1575729a2"
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[4.47696617 2.7519661  2.        ]\n"
     ]
    }
   ],
   "source": [
    "print(svd.singular_values_)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "4CbG9kJXictK"
   },
   "source": [
    "#### 非负矩阵分解"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "colab": {},
    "colab_type": "code",
    "id": "KcA2Rd4Df_DE"
   },
   "outputs": [],
   "source": [
    "def inverse_transform(W, H):\n",
    "    # 重构\n",
    "    return W.dot(H)\n",
    "\n",
    "def loss(X, X_):\n",
    "    #计算重构误差\n",
    "    return ((X - X_) * (X - X_)).sum()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "colab": {},
    "colab_type": "code",
    "id": "yRXZt6CfYPJq"
   },
   "outputs": [],
   "source": [
    "# 算法 17.1\n",
    "\n",
    "class MyNMF:\n",
    "    def fit(self, X, k, t):\n",
    "        m, n = X.shape\n",
    "        \n",
    "        W = np.random.rand(m, k)\n",
    "        W = W/W.sum(axis=0)\n",
    "        \n",
    "        H = np.random.rand(k, n)\n",
    "        \n",
    "        i = 1\n",
    "        while i < t:\n",
    "            \n",
    "            W = W * X.dot(H.T) / W.dot(H).dot(H.T)\n",
    "            \n",
    "            H = H * (W.T).dot(X) / (W.T).dot(W).dot(H)\n",
    "            \n",
    "            i += 1\n",
    "            \n",
    "        return W, H"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "colab": {},
    "colab_type": "code",
    "id": "zc1IFBBIajXk"
   },
   "outputs": [],
   "source": [
    "model = MyNMF()\n",
    "W, H = model.fit(X, 3, 200)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 125
    },
    "colab_type": "code",
    "id": "EYo7JyXXbNpJ",
    "outputId": "716dbb03-5229-4c93-a75b-81d07ff4be85"
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[0.00000000e+00, 4.27327705e-01, 6.30117924e-27],\n",
       "       [5.11680721e-97, 8.57828102e-01, 0.00000000e+00],\n",
       "       [2.97520805e-88, 2.39454414e-18, 4.36332453e-01],\n",
       "       [2.15653741e+00, 3.38756557e-21, 8.38350315e-01],\n",
       "       [7.36106818e-01, 1.00339294e-54, 8.72892573e-38],\n",
       "       [6.78344810e-01, 1.07009504e+00, 8.57259947e-01]])"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "W"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 71
    },
    "colab_type": "code",
    "id": "1cEFDsgXbnXZ",
    "outputId": "e3c8eaf0-bcd8-48f5-edee-57f643b07fed"
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[7.14647214e-10, 1.01233436e-03, 3.76097657e-02, 1.35755597e+00],\n",
       "       [9.30509415e-01, 1.86788842e+00, 1.16682319e-02, 4.54479182e-03],\n",
       "       [4.95440453e-03, 6.18432747e-04, 2.28890170e+00, 8.61836630e-02]])"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "H"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 125
    },
    "colab_type": "code",
    "id": "JFqGat0JdVlL",
    "outputId": "567c82ff-997f-49e7-b829-f9d4c828f9c9"
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[3.97632453e-01, 7.98200474e-01, 4.98615876e-03, 1.94211546e-03],\n",
       "       [7.98217125e-01, 1.60232718e+00, 1.00093372e-02, 3.89865014e-03],\n",
       "       [2.16176748e-03, 2.69842277e-04, 9.98722093e-01, 3.76047290e-02],\n",
       "       [4.15352814e-03, 2.70160021e-03, 2.00000833e+00, 2.99987233e+00],\n",
       "       [5.26056687e-10, 7.45186228e-04, 2.76848049e-02, 9.99306203e-01],\n",
       "       [9.99980721e-01, 2.00003500e+00, 2.00018226e+00, 9.99636205e-01]])"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# 重构\n",
    "X_ = inverse_transform(W, H);X_"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 35
    },
    "colab_type": "code",
    "id": "FmXCjjnyfcfY",
    "outputId": "819c8029-12d3-4344-e5af-b352b51507a3"
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "4.002356735601103"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# 重构误差\n",
    "\n",
    "loss(X, X_)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "dGu-tGDxhEcf"
   },
   "source": [
    "### 使用 sklearn 计算"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "colab": {},
    "colab_type": "code",
    "id": "sLN4FLmvb6tt"
   },
   "outputs": [],
   "source": [
    "from sklearn.decomposition import NMF\n",
    "model = NMF(n_components=3, init='random', max_iter=200, random_state=0)\n",
    "W = model.fit_transform(X)\n",
    "H = model.components_"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 125
    },
    "colab_type": "code",
    "id": "Fm5W6xQ0b_jl",
    "outputId": "d2ae79f6-0fea-47ba-a2ba-1742b78f71b1"
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[0.        , 0.53849498, 0.        ],\n",
       "       [0.        , 1.07698996, 0.        ],\n",
       "       [0.69891361, 0.        , 0.        ],\n",
       "       [1.39782972, 0.        , 1.97173859],\n",
       "       [0.        , 0.        , 0.65783848],\n",
       "       [1.39783002, 1.34623756, 0.65573258]])"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "W"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 71
    },
    "colab_type": "code",
    "id": "wheFi8ZwcCY9",
    "outputId": "9dd57216-889a-4cc4-ccb8-82e790595d59"
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[0.00000000e+00, 0.00000000e+00, 1.43078959e+00, 1.71761682e-03],\n",
       "       [7.42810976e-01, 1.48562195e+00, 0.00000000e+00, 3.30264644e-04],\n",
       "       [0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 1.52030365e+00]])"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "H"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 125
    },
    "colab_type": "code",
    "id": "9hfj3bRXgHRb",
    "outputId": "3be1a4d8-a160-4200-d4e7-f5c92f2aa1af"
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[3.99999983e-01, 7.99999966e-01, 0.00000000e+00, 1.77845853e-04],\n",
       "       [7.99999966e-01, 1.59999993e+00, 0.00000000e+00, 3.55691707e-04],\n",
       "       [0.00000000e+00, 0.00000000e+00, 9.99998311e-01, 1.20046577e-03],\n",
       "       [0.00000000e+00, 0.00000000e+00, 2.00000021e+00, 3.00004230e+00],\n",
       "       [0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 1.00011424e+00],\n",
       "       [1.00000003e+00, 2.00000007e+00, 2.00000064e+00, 9.99758185e-01]])"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "X__ = inverse_transform(W, H);X__"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 35
    },
    "colab_type": "code",
    "id": "-iALOTKfgKzP",
    "outputId": "80c73327-dd36-403f-aa20-d4e2d3f796a7"
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "4.000001672582457"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "loss(X, X__)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "----\n",
    "本章代码来源：https://github.com/hktxt/Learn-Statistical-Learning-Method\n",
    "\n",
    "本文代码更新地址：https://github.com/fengdu78/lihang-code\n",
    "\n",
    "中文注释制作：机器学习初学者公众号：ID:ai-start-com\n",
    "\n",
    "配置环境：python 3.5+\n",
    "\n",
    "代码全部测试通过。\n",
    "![gongzhong](../gongzhong.jpg)"
   ]
  }
 ],
 "metadata": {
  "colab": {
   "collapsed_sections": [],
   "name": "LSA.ipynb",
   "provenance": [],
   "version": "0.3.2"
  },
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}
